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ABSTRACT 

This  thesis  presents  a  computer  simulation  program  designed  to  predict 
the  performance  of  a  laser  designator  in  an  operational  environment.  It 
uses  aspects  of  Fourier  optical  theory,  atmospheric  modulation  transfer 
function,  atmospheric  transmittance,  and  turbulence  to  analyze  the  spot  size 
and  the  power  on  a  target.  The  guidance  system  requires  the  specific 
power  and  spot  size  on  a  target  to  track  it.  So  the  optimum  laser 
designator  for  special  military  purpose  can  be  chosen  by  the  output  of  this 
program  for  a  particular  atmospheric  condition. 
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I.      INTRODUCTION 

A.      GENERAL 

In  recent  years  as  the  naval  surface  combatants  become  more 
sophisticated,  the  use  of  laser  radiation  for  the  detection  of  targets,  and 
weapon  control  becomes  more  important  for  the  navy. 

A  laser  system  has  several  advantages  over  a  microwave  system 
because  of  its  very  short  wavelength.  The  advantages  are  higher  resolution, 
less  power  requirement  and  space  needed,  interception  and  jamming  are 
difficult,  also  it  can  be  portable.  So,  laser  guided  weapons  have  been  sold 
by  the  world's  weapon  manufacturers  to  all  nations  capable  of  supporting 
defense  forces.  But  the  laser  system  has  a  big  disadvantage  in  poor 
performance  in  bad  weather. 

The  performance  of  military  laser  systems  for  imaging,  target  detection, 
tracking,  target  designation,  warning  of  missile  launch  or  laser  irradiation, 
optical  fuzing,  or  laser  weaponry,  is  very  strongly  dependent  on  the 
transmission  (extinction)  and  modification  of  such  radiation  by  atmospheric 
effects.  The  extinction  depends  on  absorption  and  scattering  by  the 
molecular  constituents  of  the  atmosphere  and  by  the  particles  suspended  in 
it.  The  modulation  transfer  degradation,  beam  wander,  beam  broadening 
and  scintillation  effects  are  caused  by  atmospheric  turbulence.  Turbulence 
effects  over  water  as  well  as  over  land  are  a  big  problem  and  may  be  a 
limiting  factor  in  optical  systems. 


How  to  decrease  the  disadvantage  in  atmospheric  effects  and  how  to 
select  the  optimum  system  within  the  operational  environment  are  more 
important  to  a  user,  because  of  the  above  factors.  So  a  user  (such  as  the 
Korean  navy)  must  study  and  research  in  the  area  of  optical  radiation 
propagation  through  the  atmosphere.  In  support  of  this  effort,  special 
attention  to  computer  programming  of  these  effects  will  be  devoted  to  select 
an  optimum  laser  designator  within  the  operational  environment  that  includes 
the  atmospheric  transmittance  (extinction)  and  turbulence.  This  is  the 
parameter  computed  for  us  by  the  program.  The  source  of  this  algorithm 
is  based  on  the  program  written  by  Prof.   Milne  in   [Ref.   3]. 

B.      BACKGROUND 

For  a  laser  designator,  the  predicted  beam  spot  size  on  the  target  can 
be  calculated  from  the  atmospheric  modulation  transfer  function  (AMTF)  and 
line  spread  function  (LSF)  of  the  laser  output  pattern  and  system  optics,  as 
described  below. 

The  one— dimensional  (single  variable)  LSF  can  be  generated  from  the 
two— dimensional  point  spread  function  (PSF)  by  integrating  over  one  of  the 
two  dimensions  of  a  point  source  image.  This  one— dimensional  LSF  is 
easier  to  calculate  the  convolution  theorem  than  the  two— dimensional 
function  in  computer. 

The  modulation  transfer  function  (MTF)  of  the  entire  system  is  the 
AMTF  times  the  Fourier  transforms  of  the  LSF  of  the  laser  output  pattern 
and  of  the  optics  diffraction  (also  defocus  and  aberration  broadening).  This 
yields   the  convolution  of  AMTF,   source,   and  optics.      The  LSF  of  the  laser 


output  pattern  is  obtained  from  its  radial  distribution  by  integrating  over 
distance  y.     Then  the  Fourier  transform  of  the  LSF  is  calculated. 

Taking  the  inverse  Fourier  transform  of  the  above  yields  the  LSF  for 
the  entire  system.  This  LSF  can  be  converted  to  a  radial  distribution 
function  by  application  of  the  ABEL  transform.  The  fractional  energy 
falling  within  a  given  radius  is  then  calculated  by  integrating  with  respect 
to  r. 

The  algorithm,  introduced  above,  can  predict  the  size  of  a  laser  spot 
on  a  target  and  the  radient  energy  in  the  spot  for  different  systems.  So, 
we  can  compare  the  systems  to  determine  which  one  is  optimum  for  our 
operational  environment. 


II.     THEORY 

The  first  element  in  the  system  we  axe  analyzing  is,  of  course,  the 
laser  itself.  The  radial  intensity  output  of  the  laser  is  generally  Gaussan 
and  can  be  expressed  by: 

A  =   A0  exp  (-r2/2<r2) 

where         A;      maximum       intensity       of       the       laser       source       in 
watts/steradian 
r;       angular  distance  from  center  in  //rad 
<j;       angular  divergence  of  the  laser  source  in  //rad 

This  can  be  calculated  for  each  laser  for  approximate  results  and  also  can 
be  measured  directly.  One  dimensional  line  spread  function  for  the  source 
is  obtained  from  the  radial  intensity  output  pattern  of  the  laser  by 
integrating  to  obtain  the  intensity  falling  in  a  vertical  slit  of  width  delta  x 
(dx)  and  a  distance  x  from  the  center  of  the  pattern  as  shown  in  Fig  1 
[Ref  1].  A  Fourier  transform  is  then  used  on  this  LSF  to  yield  the  MTF 
of  the  laser. 

The  second  element  is  the  optics  function.  The  optics  associated  with 
the  laser  is  basically  near  diffraction  limited.  The  Airy  function  is  used  to 
calculate    the    diffraction    point    spread    function    for    the   laser   optics.      From 


[Ref    6],     the    Airy    function    for    circular    aperture    diffraction     pattern     is 
expressed  by  the  next  equation,  and  illustrated  in  Fig  2. 


IW  =  i(o) 


2  Jt(k  a  sing) 
k  a  sing 


where,        Ji(u)  =  first  order  Bessel  function  of  the  first  kind 

k        =  2tt/A 
a        =  radius  of  the  aperture 
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(a)  Radial  distribution  of  laser 
output  pattern  E(r); 
intensity     r;   radius 

(b)  Integration  to  produce  the  line 
spread  function  intensity  vs.  x 
is  the  area  of  the  hatched 
strip 

(c)  Line  spread  function  intensity 
vs.  x 


[Fig.    1]      Change  from  PSF  to  LSF 
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[Fig.   2]      Airy  function  for  circular  aperture  diffraction  pattern 

This  point  spread  function  is  then  converted  to  a  line  spread  function, 
Fourier  transformed  to  get  the  MTF  and  then  combined  by  the  convolution 
theorem  into  MTF  of  the  laser  plus  optics. 

The  third  element  to  be  calculated  is  the  effects  of  atmospheric 
turbulance  on  the  laser  beam.  A  properly  path— position  weighted  mean 
value  of  Cn  becomes  the  central  parameter  of  interest  in  calculating  the 
AMTF  due  to  turbulence.  Cn  and  Ct  are  called  respectively  the  refractive 
index  structure  constant  and  the  temperature  structure  constant,  and  are 
related  through  the  equation  below. [Ref  1] 
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Cn2    =    (        ',j,/        10-6    )2    Ct2 


Ct2   •    r2/3  =   <   (Trl  -  Tr2)2  > 


where     Tri— Tr2  ;   is   the   temperature   difference   between   two   points 

seperated  by  a  distance  r 
P       ;  pressure  in  millibars 

But  in  the  field,  the  Cn2  can  be  measured  directly  by  optical  methods 
decribed  in  [Ref  1].  The  Cn2  is  used  to  predict  the  MTF  of  the 
atmosphere  over  the  path  between  the  laser  system  and  the  target,  the 
AMTF  is  detailed  in  subsection  3. 

The  AMTF  is  then  multiplied  by  the  transmittance  (T  =  e'_//R)  ,  fi 
=  extinction  coefficient   ,  R  =  range  ) 

Now  the  total  effect  of  the  atmosphere,  laser  itself  and  optics  can  be 
found  by  the  product  of  the  three  transform  functions. (convolution  theorem) 
The  inverse  transform  then  will  yield  the  line  spread  function  of  the  total 
system,  the  LSF  as  predicted  on  the  target.  [Ref  1]  describes  the  process 
of  converting  the  one  dimensional  LSF  to  a  circular  symmetric  point  spread 
function   (PSF)  using  the  ABEL  transform,  detailed  in  subsection  2. 

After  obtaining  the  radial  intensity  distribution,  the  size  of  a  laser  spot 
on  a  target  can  be  predicted.  Also  the  power  inside  a  circle  of  radius  r 
can  then  be  obtained  by  integrating  the  PSF  out  to  r. 

A.      FOURIER  TRANSFORM  AND  CONVOLUTION 

In  order  to  get  the  total  effect  of  the  atmosphere,  laser  itself  and 
optics,  some  elements  of  Fourier  theory  must  be  considered.  The  Fourier 
transform   represents   the   one   dimensional   position   variable   g(y)    expressed   in 


the  spatial  frequency  domain  u(v).  The  inverse  transform  repeats  the 
operation  in  the  opposite  direction.  The  forms  are  commonly  represented  as 
shown. 


roo 


u(v)   =         g(y)  exp(-2?r  ivy)  dy  =  F[g(y)] 


—oo 


g(y)   =         u(v)  exp(27r  ivy)  dv  =  F  1[u{v)} 

The  program  calculates  the  value  of  the  integral  at  a  preset  number  of 
points  and  yields  a  discrete  Fourier  transform.  A  theorem  from  Fourier 
theory  that  allows  us  to  compute  the  line  spread  function  for  one  of  the 
elements  in  the  system  given  is  the  convolution  theorem  stated  as  follows. 

If  F[g(y)]   =   G(v),  and 

F[h(y)]   =   H(v) 
then  F[g*h]   =  GH         (*  represents  the  convolution) 

B.      ABEL  TRANSFORM   [Ref.l] 

In  order  to  get  circularly  symmetric  two— dimensional  PSF  from 
one— dimensional  LSF  the  ABEL  transform  [Griem,  1964;  Magnus,  1949]  is 
used  in  this  program  referred  to  above  in  section  II.  The  one  dimensional 
line  spread  function  can  be  utilized  to  calculate  the  corresponding  radial 
distribution  of  intensity  in  the  two— dimensional  image,  provided  that  the 
intensity  is  known  to  be  circularly  symmetric,  i.e.  the  intensity  E(r)  in  the 
image    depends    only   on    the   radial    distance   from    the   center   of   the    image. 
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This  inversion  is  called  the  ABEL  transform.  For  such  an  image,  the  line 
spread  function  (LSF)  is  obtained  from  the  radial  distribution  by  integrating 
the  intensity  lying  within  an  infinitesmal  slit  of  width  dx. 


I(x)  =     j°°     E(r)  dy 


=   2 


00 


r-R 


0 


E(  J    x2  +  y2     )  dy 


where     R  is  a  radius  beyond  which  there  is  no  intensity 

The   ABEL   transform  allows   us   to   reverse  this   process.      It   is   accomplished 
by  earring  out  the  calculation  indicated  by; 


E(r)   =  - 


R  T/ 


r(x)  dx 

2_v2^l/2 


r      (r2-x2) 


where  I(x)   =   0  for  x   >   R 

I'(x)  is  the  x  derivative  of  I(x) 

This   integral   can   be   carried   out   by   step   wise   numerical   integration   with   a 
computer. 

C.      MTF  OF  ATMOSPHERE 

The  modulation   transfer  function  is  a  function  of  the  spatial  frequency, 
angular     spatial     frequency     f     given     in     (cycles/radian)     or     linear     spatial 


frequency  v  given  in  (cycles/meter).  The  one— dimensional  MTF  is  directly 
applicable.  Images  are,  of  course,  really  two— dimensional,  and  much  work 
has  been  done  in  the  field  of  two— dimensional  image  processing,  both  analog 
and  digital.  However,  for  many  purposes,  the  one— dimensional  form  is 
adequate,  and  is  much  more  readily  handled  analytically  and  in  computer 
processing.  We  already  saw  that  for  circularly  symmetric  images,  the 
two— dimensional  behavior  can  be  obtained  from  the  one— dimensional.  The 
modulation  transfer  function  has  the  added  advantage  that  for  a  total 
optical  system  the  system  MTF  is  the  product  of  the  MTF's  of  the 
components.  So,  the  MTF  of  atmosphere  can  be  expressed  independently  by 
its  own  parameters. 

There  are  two  limiting  situations  used  in  calculating  the  AMTF. 
First,  the  short  exposure,  or  short  term  MTF  is  obtained  by  taking  the 
average  of  successive  LSF  in  which  each  LSF  is  shifted  such  that  the  center 
of  the  area  of  each  curve  has  a  common  origin.  Also  it  discribes  images 
obtained  in  the  limit  of  time  intervals  sufficiently  short  that  the  turbulence 
can  be  considered  frozen.  The  average  short  term  MTF  is  then  the  average 
of  the  MTF's  corresponding  to  many  such  images.  This  procedure,  in 
effect,  removes  the  image  wander.  Second,  the  long  exposure,  or  long  term 
MTF  is  obtained  by  including  the  image  wander,  by  not  performing  a  shift 
in  the  origin.  Also  it  results  simply  from  taking  an  integrated  image  of  a 
sufficiently  long  time  which  sees  effectively  all  possible  turbulence 
configurations.  It  is  seen  to  be  the  result  of  image  spread  (the  short  term 
MTF)  plus  image  wander.  So  the  short  term  LSF  is  narrower  than  the 
long  term  LSF. 
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From     the     [Ref     1]     AMTF     can     be     summarized     by     the     following 
expression. 


M(f)  =  exp 


-57.64  0  Cn2  Z0  A'1/3  f5/3  [l-a  (-^-)l/3] 


where         A        =  wavelength 

f        =  angular  spatial  frequency(cycles/radian) 

D       =  diameter  of  optics 

Zo      =  range 

a       =  0  for  the  L.T.   MTF 

=  1/2  for  the  far-field  S.T.  MTF,  D<<(ZA)1/2 

=  1  for  the  near-field  S.T.  MTF,  D^ZA)1/2 

0       =  1   for  a  plane  wave 

=  3/8  for  a  spherical  wave 
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III.     PROGRAM  DESCRIPTION 

This  program  can  be  run  on  IBM  personal  computers  or  compatibles 
under  the  DOS  operating  system  (version  2.0  or  later),  and  consists  of 
several  subroutines,  a  function  and  a  mainbody.  All  program  modules  are 
coded  in  MS's  ANSI  77  Fortran  Compiler.  The  source  file  is  compiled  and 
linked  to  library  with  command  (FL  /G2),  then  the  .EXE  file  is  obtained 
and  this  program  is  ready  to  be  run. 

The  source  file  is  listed  in  Appendix  A,  the  typical  input  data  file  and 
module  interconnections  are  listed  in  Appendix  B. 

Some  input  data  and  modules  are  added  to  the  program  written  by 
Prof.  Milne  in  [Ref.  3]  for  some  parameters  to  be  varied,  and  for  the 
desired  output  to  be  taken. 

A.      INPUT  DATA 

The  general  input  data  to  this  program  must  be  entered  interactively 
from  the  keyboard  or  from  the  input  file.  In  either  case,  these  data  must 
be  determined  by  the  user  prior  to  running  this  program.  The  input  data 
are  described  and  determined  briefly  below. 

1.       ID  ATE 

Time  and  date  of  the  run,  60  characters,  it  will  be  placed  in  the 
header  of  the  output  file. 
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2.  INPT 

This  input  indicates  the  type  of  laser  source  function  to  be  used 
in  the  program.  If  1  is  entered,  then  a  gaussian  source  function  is 
calculated.  If  2  is  entered,  then  the  source  function  is  inputted  from  the 
input  file.  If  3  is  entered,  then  combination  function  of  source-optics  is 
entered  from  the  input  file. 

3.  IHEAD(l) 

Information  about  the  laser  source,  80  characters,  it  will  be  placed 
in  the  header  of  the  output  file. 

4.  WAV 

Wavelength  of  the  laser  in  meters,  Format  E10.4,  this  input  can 
come  from  the  characteristic  of  the  laser  designator  determined  by  the 
weapon  manufacturer. 

5.  OBRA 

Obscuration  ratio,  Format  F7.5,  the  ratio  of  the  diameter  of  the 
inner  mirror  to  outer  mirror.      This  data  can  come  from  manufacturer  also. 

6.  RANGE 

The  distance  between  the  laser  designator  and  the  target  in 
meters,  Format   F7.0. 

7.  OPT 

Diameter  of  the  optics  in  meters,  Format   F7.3. 

8.  SQCN 

Squared  refractive  index  structure  constant  in  meters^"2'  3',  Format 
E10.4.  This  input  will  be  the  path— weighted  integrated  value  in  the 
operational  environment  and  can  be  measured  directly  by  optical  methods 
[Ref.   1]. 
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9.  SCALE 

Angular  scale  in  microradians  per  point,  Format  F7.0. 

10.  EXT 

The  extinction  coefficient  of  the  operational  atmosphere  in 
meters"1,  Format  E10.4,  this  data  can  come  from  Lowtran  6  in  NPS  IBM 
3033AP  computer. 

11.  INOP 

Input  type  of  the  optics  diffraction  pattern  to  be  used  in  the 
program.  If  1  is  entered,  then  the  optics  diffraction  function  is  calculated. 
If  2  is  entered,  then  it  is  entered  from  the  input  file. 

12.  IHEAD(2) 

Information  about  the  optics,  80  characters,  it  will  be  placed  in 
the  header  of  the  output  file. 

13.  A 

The  maximum  intensity  of  the  laser  source  in  watts/steradian, 
Format   E10.4,  this  data  can  come  from  manufacturer  also. 

14.  B 

Sigma(standard  deviation)  of  the  laser  source  (gaussian)  in 
microradians,  Format  El 0.4.  This  data  can  be  determined  to  measure  the 
angle  between  the  beam  center  and  the  line  at  which  the  laser's  intensity 
falls  to  e-1'2  of  its  peak  intensity.  This  data  can  come  from  manufacturer 
also. 

15.  EXPO 

The  type  of  exposure  pattern  to  be  taken  for  output.  If  1  is 
entered,  then  the  long— term  (exposure)  is  taken.  If  2  is  entered,  then  the 
short— term   (exposure)  is  taken. 
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16.  FRAP 

Fraction  of  total  power  in  watts  that  is  required  at  the  target. 

17.  VARY 

Selection  of  the  run.  If  1  is  entered,  then  one  of  input 
parameters  will  be  varied  in  program.  If  2  is  entered,  then  fixed  value  of 
input  parameters  will  be  used  in  program. 

18.  SEL 

Selection  of  the  parameter  that  will  be  varied  in  program. 
SQCN,  EXT,  RANGE,  and  OPT  can  be  chosen  to  be  varied.  The 
maximum  and  minimum  of  that  value  should  be  entered.  Input  Format  is 
same  as  above  6,   7,  8,  and   10. 

B.      MODULE  DESCRIPTIONS 
1.       GETVAR 

a.  Module  type  :      Fortran  77  Subroutine 

b.  Usage  :     Call  GETVAR(min,max,va) 

c.  Function  :  To  get  20  points  between  minimum  and 
maximum 

d.  Interface  : 

(1)  External  references  :     Float 

(2)  Common  data  :     None 

(3)  Formal  arguments   : 

(a)  Min  (real),  minimum  of  the  variable. 

(b)  Max  (real),  maximum  of  the  variable. 

(c)  Va  (real   array),   20   values  of  the  variable  between 
minimum  and  maximum. 
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(4)  Input   :     Min,  Max 

(5)  Output   :     Va 

2.  AIRY 

a.  Module  type  :     Fortran  77  Real  function 

b.  Usage  :     Y  =  AIRY  (X) 

c.  Function   :      Calculates   the  value  of  the  Airy  function  of  x. 
AIRY(x)   =  2  *  Ji(x)/x,  where  Ji(x)  is  the  Bessel  function  of  order  one. 

d.  Interface  : 

(1)  External  references  :     Float,  ABS,  MAX1 

(2)  Common  data  :     None 

(3)  Formal  arguments  : 

(a)     X   (real),   value  of  the   independent   variable  of  the 
Airy  function 

(4)  Input   :      X 

(5)  Output   :      Airy 

3.  LSF 

a.  Module  type  :     Fortran  77  Subroutine 

b.  Usage  :     Call  LSF(RG3,k) 

c.  Function    :       Converts    a    point    spread    function    into    a    line 
spread  function. 

d.  Interface  : 

(1)  External  references  :     Float,  SQRT,   IFIX 

(2)  Common  data  :      None 

(3)  Formal  arguments  : 

(a)     RG3    (real    array),    this    is    the    PSF    on    input    and 
the  LSF  on  output. 
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(b)     K   (integer),  number  of  points  in  the  array  RG3. 

(4)  Input   :      RG3,  K 

(5)  Output   :      RG3 

4.  AMTF 

a.  Module  type  :     Fortran  77  Subroutine 

b.  Usage  :     Call  AMTF(j) 

c.  Function  :  Calculates  the  atmospheric  MTF  from  given 
value  for  long— term  and  short— term,  and  then  the  result  is  multiplied  by 
extinction  factor  (e~^   ). 

d.  Interface  : 

(1)  External  references  :     EXP,  COS 

(2)  Common  data   :      /real/   RGl,   RG2.      /   /   WAV,   OPT, 
OBRA,   SCALE,   RANGE,  SQCN,  EXT,  A,  B,  IDATE,   IHEAD 

(3)  Formal  arguments  : 

(a)  J  (integer),  if  zero  on  input,  then  the  long  term 
AMTF  is  calculated;  if  one  on  input,  then  the  short  term  AMTF  is 
calculated. 

(4)  Input    :        J,    SQCN,    RANGE,    WAV,    EXT,     SCALE, 
OPT 

(5)  Output   :      RG2 

5.  FFT 

a.  Module  type  :     Fortran  77  Subroutine 

b.  Usage  :     Call  FFT(RG,ID,KA) 

c.  Function  :  Calculates  the  Fourier  transform  or  inverse 
Fourier  transform  of  a  real— value  symmetric  function  using  the  Fast  Fourier 
transform  (FFT)  algorithm. 
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d.       Interface  : 

(1)  External  references  :  Float,  CMPLX,  CEXP,  CONJG 

(2)  Common  data  :     None 

(3)  Formal  arguments  : 

(a)  RG  (real),  on  input,  RG  is  an  array  of  functional 
values;  on  output,  it  is  the  forward  (or  inverse)  Fourier  transform  of  these 
functional  values. 

(b)  ID  (integer),  specifies  whether  forward  (ID  =  1)  or 
inverse  (ID  =  0)   Fourier  transform  is  to  be  calculated. 

(c)  KA  (integer),  the  largest  non-zero  term  of  forward 
(or  inverse)   Fourier  transform. 

(4)  Input   :      RG,   ID 

(5)  Output   :     RG,  KA 
6.       BITINV 

a.  Module  type  :     Fortran  77  Subroutine 

b.  Usage  :      Call  BITINV(K,L,M) 

c.  Function   :      Finds   the  bit   inversion   of  an   M   bit   number   K 
to  L.      This  is  used  in  the  Fourier  transform  subroutine. 

d.  Interface  : 

(1)  External  references  :  None 

(2)  Common  data  :      None 

(3)  Formal  arguments  : 

(a)     K   (integer),   the  integer  whose  bit   inverse  is   to  be 


found. 


(b)     L  (integer),  the  bit  inverse  of  K. 
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(c)     M   (integer),   the  number  of  bits   in   the  integers   K 
and  L. 

(4)  Input   :      K,   M 

(5)  Output   :     L 

7.  ABEL 

a.  Module  type  :     Fortran  77  Subroutine 

b.  Usage  :     Call  ABEL(K) 

c.  Function  :  Converts  a  line  spread  function  (x)  to  a  point 
spread  function  (r). 

d.  Interface  : 

(1)  External  references  :  Float,  SQRT,  ABS 

(2)  Common  data  :      /real/  RG1,   RG2 

(3)  Formal  arguments  : 

(a)     K   (integer),   number  of  points   in   array   RG2   to  be 
processed. 

(4)  Input   :      RG2,   K 

(5)  Output   :      RG2 

8.  INT 

a.  Module  type  :      Fortran   77  Subroutine 

b.  Usage  :     Call  INT(K,TPOW) 

c.  Function  :  Integrates  the  point  spread  function  to  find  the 
total  power  as  a  function  of  radius  reached  around  the  target  and 
normalizes  it. 

d.  Interface  : 

(1)     External  references  :   Float 
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(2)  Common  data  :      /real/   RG1,   RG2.      /   /   WAV,   OPT, 
OBRA,   SCALE,  RANGE,  SQCN,  EXT,  A,  B,  IDATE,  IHEAD 

(3)  Formal  arguments  : 

(a)  K  (integer),  number  of  points  in  array  RG2. 

(b)  TPOW  (real),  total  integrated  power  reached 
around  the  target. 

(4)  Input   :     K,  RG2,  SCALE,  RANGE 

(5)  Output   :      TPOW 
9.       OUTIN 

a.  Module  type  :      Fortran  77  Subroutine 

b.  Usage  :     Call  OUTIN(BUF,J,K,L) 

c.  Function   :      Control   the  input   and  output   files   to   read   and 
write  the  data. 

d.  Interface  : 

(1)  External  references  :     None 

(2)  Common    data    :       /OI/    LU,    LI,    IH,    IM,    IS,    IHS,    IY, 
IMO,  ID 

(3)  Formal  arguments  : 

(a)  BUF  (array),  array  data  which  are  to  be  read  or 
written. 

(b)  J  (integer),  if  -1  on  input  then  call  the  OPEN 
subroutine;  if  0  on  input  then  call  the  READF  subroutine;  if  1  on  input 
then  call  the  WRITF  subroutine. 

(c)  K  (integer),  the  count  number  of  variable  to  be 
varied. 
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(d)     L    (integer),   parameter   specifying    which   title   is    to 
be  printed  on  output  file  or  which  data  is  to  be  read  from  input  file. 

(4)  Input   :      BUF,  J,  K,  L 

(5)  Output   :     BUF 

10.  OPEN 

a.  Module  type  :      Fortran  77  Subroutine 

b.  Usage  :     Call  OPEN 

c.  Function  :     Open  the  input  and  output  files. 

d.  Interface  : 

(1)  External  references  :     None 

(2)  Common   data    :       /OI/    LU,    LI,    IH,    IM,    IS,    IHS,    IY, 
IMO,   ID 

(3)  Formal  arguments  :      None 

(4)  Input   :      None 

(5)  Output   :      None 

11.  READF 

a.  Module  type  :     Fortran  77  Subroutine 

b.  Usage  :     Call  READF(IER,BUF,L) 

c.  Function  :      Read  the  data  from  input  file. 

d.  Interface  : 

(1)  External  references  :     None 

(2)  Common   data   :      /OI/    LU,    LI,    IH,    IM,    IS,    IHS,    IY, 
IMO,  ID 

(3)  Formal  arguments  : 

(a)     IER    (integer),    error    number    of   reading    the    input 
data  file. 
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(b)  BUF  (array),  array  data  which  are  to  be  read 
from  input  data  file. 

(c)  L  (integer),  if  1  on  input  then  read  the 
information;  if  2  on  input  then  read  the  source  function  or  source— optics 
function. 

(4)  Input   :      L 

(5)  Output   :      IER,  BUF 

12.  WRITF 

a.  Module  type  :     Fortran  77  Subroutine 

b.  Usage  :     Call  WRITF(IER,BUF,L) 

c.  Function  :     Write  the  data  in  output  file. 

d.  Interface  : 

(1)  External  references  :     None 

(2)  Common   data    :       /OI/    LU,    LI,    IH,    IM,    IS,    IHS,    IY, 
IMO,   ID 

(3)  Formal  arguments  : 

(a)  IER  (integer),  error  number  of  writing  the  output 
data  file. 

(b)  BUF  (array),  array  data  which  are  to  be  written 
on  output  data  file. 

(c)  L  (integer),  parameter  specifying  which  title  is  to 
be  printed  on  output  file. 

(4)  Input   :      BUF,   L 

(5)  Output   :      IER 

13.  RESUL 

a.       Module  type  :      Fortran  77  Subroutine 
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b.  Usage  :     Call  RESUL(INU,TPOW,RADIU,POWER,FRAP) 

c.  Function    :      Calculates    the   radii   on   the   target   and   average 
usable  power. 

d.  Interface  : 

(1)  External  references  :  Float 

(2)  Common  data  :      /real/   RGl,   RG2.      /   /   WAV,   OPT, 
OBRA,   SCALE,  RANGE,  SQCN,  EXT,  A,  B,  IDATE,   IHEAD 

(3)  Formal  arguments  : 

(a)  INU   (integer),   the  count   number  of  variable  to  be 
varied. 

(b)  TPOW     (real),     total     integrated     power     reached 
around  the  target. 

(c)  RADIU    (real    array),    radii    on    the    target    which 
depends  on  the  variable  VA. 

(d)  POWER   (real    array),    average   usable   power    which 
depends  on  the  variable  VA. 

(e)  FRAP     (real),     function     of    total     power     that     is 
required  at  the  target. 

(4)  Input   :      RG2,  FRAP,  RANGE,  SCALE,  TPOW 

(5)  Output   :      POWER,  RADIU 

14.     WRITO 

a.  Module  type  :     Fortran  77  Subroutine 

b.  Usage  :     Call  WRITO(VA,RADIU,POWER,ABC) 

c.  Function  :     Write  the  data  in  OUT2.  file. 

d.  Interface  : 
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(1)  External  references  :  None 

(2)  Common  data  :     /OJ/  SEL,  VERY 

(3)  Formal  arguments  : 

(a)  VA  (real  array),  20  values  of  the  variable  between 
minimum  and  maximum. 

(b)  RADIU  (real  array),  radii  on  the  target  which 
depends  on  the  variable  VA. 

(c)  POWER  (real  array),  average  usable  power  which 
depends  on  the  variable  VA. 

(d)  ABC  (array),  information  that  includes  the 
common  data  /  /  WAV,  OPT,  OBRA,  SCALE,  RANGE,  SQCN,  EXT,  A, 
B,   IDATE,  IHEAD. 

(4)  Input   :     VA,  RADIU,  POWER,  ABC,SEL 

(5)  Output   :     None 
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IV.      USE  OF  PROGRAM 

The  purpose  of  this  simulation  is  to  use  the  program  to  predict  the 
atmospheric  effects  on  a  laser  designator,  and  to  select  the  optimum  laser 
designator  within  the  operational  environment.  This  simulation  program 
includes  only  one  leg  from  a  laser  designator  to  a  target,  but  the  total 
infrared  system  including  laser  designator  must  be  described  before  simulating 
with  this  program. 

The  total  infrared  system  including  laser  designator  with  definite 
military  mission  can  be  represented  by  the  following  block  diagram  shown  in 
Fig.   3.     The  specific  function  of  the  blocks  are  discussed  in   [Ref.   7]. 

This  total  system  involves  double  pass  operation  from  a  transmitter 
(laser  designator)  to  a  target  and  back  to  the  receiver.  The  received  signal 
may  be  used  to  determine  the  position,  range,  or  velocity  of  the  target, 
then  the  guidance  system  can  be  directed  at  the  target.  A  target 
designator  system  is  the  first  leg  of  this  total  system  in  which  the  receiver 
is  physically  removed  from  the  transmitter.  The  most  obviously  important 
measure  of  system  performance  is  the  maximum  range  of  operation,  or  the 
ratio  of  received  power  to  transmitted  power.  The  range  equation  for  this 
double  pass  is  stated  by  following  expression  in  [Ref.   7]. 


p     -        Pt  a   Vt  A0  A  t 
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where         Pt      =       output    power   of   the   transmitter    (laser    designator) 

in  watts. 
Pr      =       received  power  of  the  detector  in  watts. 
At      =       aperture  of  the  transmitter  in  meters. 
A0     =       aperture  of  the  detector  in  meters. 
7/t       =       overall    path   efficiency    (=    Td    *    Tm    *    T0,    Td   is 
atmospheric      transmittance      from      designator      to 
target,   Tm  is  atmospheric  transmittance  from  target 
to  detector,  T0  is  optical  transmittance). 
Rd     =       distance  of  designator  from  target. 
Rm     =       distance  of  detector  from  target. 
a       =       optical     reflection     efficiency     of     the     target     (the 
product     of    spot     area    on     target     and     reflection 
efficiency  of  target  surface  and  cos#). 
This    equation    based   on   some   assumptions,    the   designator    (transmitter)    and 
detector  are  on   different   path,   so  the  angle  between  designator  and  detector 
on   target   is    0,   the  target   acts   as   a   diffuse   radiator,   and   the   spot   size   on 
the  target  is  smaller  than  the  target  size. 

From  the  assumptions,  the  reflected  flux  per  unit  solid  angle  for  a 
perfectly  diffuse  reflector  is  proportional  to  the  cosine  of  the  angle  between 
direction  of  interest  and  the  normal  to  the  surface,  a  relationship  known  as 
Lambert's  cosine  law.[Ref.  8] 

For  the  guidance  system  shown  in  lower  part  of  Fig.  3,  there  exists  a 
minimum  IR  power  input  from  the  target  necessary  for  detection.  This 
must  be  in  the  form  of  a  target  irradiance  at  the  optical  receiver  in 
guidance  system.      The  minimum  detectable  spectral  intensity  from  the  target 
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depends  on  the  distance  between  the  target  and  the  optical  receiver,  and 
minimum  spectral  power  required  in  the  receiver  that  includes  the  optical 
transmission  factor  and  spectral  responsivity  of  detector.  So  the  range 
equation  for  guidance  system  can  be  constructed  to  get  the  minimum 
detectable  spectral  intensity  from  the  target  by: 


t    t  \  \  *r( A)    *   R 

Irt^       ;        -MA)-R 


A0     •   e 


where 


Irt(A) 

Pr(A) 

A0 

MA) 

R 


spectral       intensity       reflected      from       target 

(Watts/sr) 

minimum    spectral    power    required    in    receiver 

(Watts) 

aperture  of  the  detector  (m2) 

spectral  extinction  coefficient   (m_1) 

range  between  target  and  guidance  system 


Finally  the  required  spectral  intensity  on  target  from  the  laser 
designator  can  be  achieved  considering  the  optical  reflection  efficiency  of 
target. 


Iit(A)  =  -^L 


where         Iit(A) 


required   spectral   intensity   on   target   from   the 

laser  designator  (Watts/sr). 

optical  reflection  efficiency  of  the  target. 
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In  order  to  track  a  target,  the  tracker  or  the  guidance  system  must  be 
directed  at  the  target  within  a  certain  angle,  and  the  signal  to  noise  ratio 
of  the  target  signal  must  be  large  enough  that  the  tracker  can  detect  the 
target.  Thus  the  prime  military  requirement  for  a  laser  designator  is 
optical  power  on  the  selected  target.  The  required  spectral  intensity  on 
target  must  be  achieved  in  a  variety  of  operating  conditions. 

And  now,  the  first  leg  of  this  total  system  shown  in  upper  part  of 
Fig.  3  can  be  treated  to  select  the  optimum  designator  in  operational 
environment,  and  was  decribed  in  Chapter  II. 


28 


0) 


Later 
Dealfaator 


l    Target    Y 


AtteaaatlBf 
Ataioapkere 


a> 


Back  f  roana* 
v     Nolae     y 


Optical 
Recetrer 


1 


Optical 

MoOlator 


I 


Detector 
(B  Co-oil  a|  Syatenl 


I 


Slfaal 
Proceaalag 


Coatrol 
Llak 


Galiaace 
System 


[Fig.   3]      Block  diagram  of  total  infrared  system 
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A.      METHOD   OF  SELECTION 

The  typical  values  of  laser  designators  and  guidance  systems  are 
classified  and  are  not  allowed  to  public  use.  For  the  purpose  of  this 
simulation,  the  assumed  data  is  to  be  used  or  only  the  method  for  test  is 
to  be  explained. 

Also  the  one  sigma  (la)  point  (40%)  of  total  integrated  power  is  to 
be  used  to  predict  the  spot  size  on  target  place.  The  value  of  sigma 
represents  the  inflection  points  of  the  Gaussian  laser  power  output  (the 
points  at  which  the  curve  changes  from  turning  downward  to  turning 
upward),  and  40%  of  the  output  power  is  included  within  these  points. 

The  averaged  usable  power  on  target  is  to  be  used  to  predict  the 
available  power  for  tracking  and  guidance  system.  This  power  is  more 
useful  considering  pointing  accuracy  and  tracking  stability. 

Four  parameters,  Cn2,  Extinction  coefficient,  Diameter  of  optics,  Range 
are  available  to  be  varied  in  this  simulation  program. 

Now  some  examples  of  using  this  method  for  this  simulation  program 
are  to  be  introduced.  First,  how  is  the  spot  size  changed  when  Cn2  is 
varied.  Second,  how  is  the  usable  power  changed  when  extinction  coefficient 
is  varied.  Third,  how  is  the  spot  size  changed  when  range  is  varied. 
Fourth,  how  can  the  system  parameter  be  selected  when  specific  operational 
environment  is  given  (single  value  of  Cn2).  Finally,  how  can  the  optimum 
system  be  selected  when  special  operational  environment  is  given  (various 
values  of  Cn2). 

1.       TURBULENCE  EFFECTS  (RUN   1) 
a.       Input  data 

(1)      Wavelength   ;    1.06  /mi 
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target 


(2)  Obscuration  ratio  ;  0.00 

(3)  Diameter  of  optics  ;  0.95xl0"i  m 

(4)  Scale  ;   1.0  ^rad/point 

(5)  Max.  intensity  ;  500  watts/steradian 

(6)  Sigma  ;   1.5  //rad. 

(7)  Cn2  ;   8.885x10-15  -  2.258xl0'i3  m^/s 

(8)  Range  ;   10  Km 

(9)  Extinction  coefficient   ;  0.5xl0*4  nr1 

(10)  Long— term  (exposure)  is  to  be  taken 

(11)  One  sigma  point  is  to  be  taken  for  fractional   power  on 


b.       Output   ;  see  Appendix  C,  A. 

2.  EXTINCTION  EFFECTS   (RUN  2) 

a.  Input  data  ;  same  as  RUN   1  except   (7)  and  (9) 
(7)     Cn2  ;   l.OxlO'14  m-2/3 

(9)     Extinction  coefficient   ;  0.5xl0"5  -  0.5x10-3  nr1 

b.  Output   ;  see  Appendix  C,   B. 

3.  RANGE  EFFECTS   (RUN  3) 

a.  Input  data  ;  same  as  RUN   1  except   (7)  and  (8) 

(7)  Cn2  ;   l.OxlO-14  m-2/3 

(8)  Range  ;  5  -  15  Km 

b.  Output   ;  see  Appendix  C,  C. 

4.  SELECTION  OF  DIAMETER  (RUN  4) 

a.       Input  data  ;  same  as  RUN   1  except   (3)  and  (7) 
(3)     Diameter  of  optics  ;  0.25X10"1  -  1.25xl0"i  m 
(7)     Cn2  ;   1.0x10-14  m-2/3 
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b.       Output  ;  see  Appendix  C,  D. 
5.       COMPARE  TWO  SYSTEMS 

Here  two  individual  system  parameters  like  (1)  through  (6)  in 
RUN  1  are  to  be  used  for  each  run.  But  the  environment  parameters  like 
(7)  through  (10)  in  RUN  1  are  to  be  taken  in  the  same  values,  because 
these  two  systems  must  be  compared  in  same  atmospheric  environment. 

Also  real  spectral  intensity  required  on  target  [Iit(A)]  instead  of 
one  sigma  point  must  be  used  to  determine  which  one  is  more  optimum  for 
the  purpose  including  tracking  or  guidance  system. 

As  a  result  the  output  will  show  spot  sizes  and  fractional  powers 
inside  a  circle  for  each  run,  then  the  user  must  make  a  decision  which 
system  is  the  best  for  their  purpose  and  operational  environment. 

B.      RESULTS  &  CONCLUSION 

This  program  can  be  useful  in  predicting  the  performance  of  laser 
designator  in  operational  environment.  It  was  used  to  show  how  the  spot 
size  and  the  averaged  usable  power  are  changed  when  some  parameters  are 
varied  as  shown  in  above  section  A. 

In  Run  1,  the  Cn2  was  varied.  The  Cn2  values  were  taken  for 
Monterey  bay  area  from  Ref.  3.  As  a  result  the  radius  of  spot  size  was 
changed  from  0.24m  to  1.36m  as  well  as  the  averaged  usable  power  was 
changed  from  .997  to  .2726xl0"1watts/sr. 

In  Run  2,  extinction  coefficient  was  varied.  The  radius  of  spot  size 
was  not  changed  but  the  averaged  usable  power  was  changed  from 
.1411xl0+1   to   .9998x10-2   watts/sr.      From  the   results   of  Run   1    and   Run   2, 
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it  is  shown  that  the  performance  of  a  laser  designator  is  affected 
significantly  by  atmospheric  effects. 

In  Run  3,  range  was  varied.  As  a  result  the  radius  of  spot  size  was 
changed  from  .09  to  .465m  as  well  as  the  averaged  usable  power  was 
changed  from  2.234  to  .4551  watts/sr.  From  the  results,  the  maximum  and 
minimum  range  can  be  predicted  in  specific  environment. 

It  was  used  to  select  the  optimum  diameter  of  optics  as  shown  in 
Run  4.  The  variation  was  not  too  big,  and  over  .83cm  of  diameter  the 
radius  of  spot  size  and  averaged  usable  power  was  not  changed.  From  the 
results,  the  optimum  diameter  of  optics  can  be  predicted. 

This  program  can  also  be  used  to  select  the  optimum  system  in  special 
operational  environment  and  purpose  like  the  Pacific  Ocean,  the  Atlantic 
Ocean,  and  air— to— air,  air— to— surface,  surface— to— surface  systems. 

In  the  military  use,  designators  are  used  to  mark  or  designate  a 
target,  thereby  turning  the  guidance  system  into  a  beacon  of  reflected 
energy  from  target.  Given  this  use,  the  required  power  on  the  target  must 
be  achieved  over  a  variety  of  operating  environments.  Also  laser  system 
user  must  look  at  such  factors  as  target  size  and  target  cross  section, 
pointing  accuracy  and  tracking  stability  requirement,  as  well  as  atmospheric 
conditions  like  turbulence  and  attenuation.  Designators  will  have  energy 
amplitude  stability,  side— lobe,  and  radiation  outside  the  main  beam 
requirements  to  avoid  loss  of  lock— on  or  misdirection  of  seekers  sensing  the 
reflected  energy.  Also  no  visible  light  can  be  emitted  and  the  beam  width 
must  be  reduced  to  minimize  countermeasure  effectiveness. 
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APPENDIX     A 


PROGRAM  LISTING 

cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 
c  c 

c  c 

C  THIS  PROGRAM  PREDICTS  THE  PERFORMANCE  C 

C  OF  THE  LASER  DESIGNATOR  WITHIN  C 

C       OPERATIONAL  ENVIRONMENT  TO  SELECT  THE  OPTIMUM  SYSTEM.        C 
C  C 

C  C 

CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 

INTEGER*2  IDATE , IHEAD , INFO , IH , IM , IS , IHS , IY , IMO , ID , LI , LU, SEL , VARY 

C0MM0N/REAL/RG1(512),RG2(512) 

C0MM0N//WAV, OPT, OBRA, SCALE, RANGE, SQCN, EXT, A, B,IDATE(30) , IHEAD (80) 

COMMON/OI/LU,LI,IH,IM,IS,IHS,IY,IMO,ID 

COMMON/OJ/SEL,VARY 

DIMENSION  INF0(128) ,ABC(64) ,VA(20) ,RADIU(20) ,P0WER(2O) 

EQUIVALENCE  (WAV, INFO (1) I, (INFO(l) ,ABC(1)) 
C  DATA  INPUT  SECTION 

LI=7 

LU=8 

CALL  getdat(iy,imo,id; 

CALL  OUTIN (ABC, -1,1,1 

CALL  OUTIN(ABC, 0,1,1)' 

WRITE(* , 1000) (IDATE(I) , 1=1 ,30) 
1000  FORMATC  ENTER  TIME  AND  DATE  OF  THE  RUN  -  60  CHARACTERS  '/1X,30A2, 
$/) 

READ(*,1030)(IDATE(I),I=1,30) 
10  WRITE(*,1010) 
1010  FORMAT(/'  ENTER  1  IF  SOURCE  IS  GAUSSIAN  TO  BE  CALCULATED '/'  ENTER 
$2  IF  THE  SOURCE  IS  TO  BE  INPUTTED '/'  ENTER  3  IF  THE  SOURCE- OPTICS 
SFUNCTION  IS  TO  BE  INPUTTED    •) 

READ(*,*)INPT 

WRITE(*,1020)(IHEAD(I),I=1,40) 
1020  FORMAT(/'  ENTER  SOURCE  INFORMATION  -  80  CHARACTERS'/ IX, 40 A2,/) 

READ(*,1030)(IHEAD(I),I=1,40) 
1030  F0RMAT(4OA2) 

WRITE(*,1040)WAV 
1040  FORMAT(/'  ENTER  WAVELENGTH  IN  METERS' ,2X,E10. 4,/) 

READ(*,*)WAV 

WRITE(*, 1050) OBRA 
1050  FORMAT(/'  ENTER  OBSCURATION  RATIO' ,2X,F7.5,/) 

READ (*,*) OBRA 

WRITE(*,1060)RANGE 
1060  FORMAT(/'  ENTER  RANGE  IN  METERS' ,2X,F7.0,/) 


34 


READ (*,*) RANGE 

WRITE(*, 1070) OPT 
1070  FORMAT(/'  ENTER  DIAMETER  OF  OPTICS  IN  METERS' ,2X,F7. 3,/) 

READ(*,*)OPT 

VRITE(*,1080)SQCN 
1080  FORMAT(/'  ENTER  CN  SQUARED  IN  METERS** (-2/3) ' ,2X,E10.4,/) 

READ(*,*)SQCN 

VRITE(*,1090)SCALE 
1090  FORMAT(/'  ENTER  SCALE  IN  MICRORADIANS  PER  POINT' ,2X,F7.2,/) 

READ(*,*)SCALE 

VRITE(*,1100)EXT 
1100  FORMAT(/'  ENTER  EXTINCTION  COEFFICIENT  IN  (1/METERS) ' ,2X,E10.4,/) 

READ(*,*)EXT 

IF(INPT.EQ.3)  GO  TO  11 

WRITE?*, 1015) 
1015  FORMATf/'  ENTER  1  IF  OPTICS  DIFFRACTION  PATTERN  IS  TO  BE  CALCULATE 
$D'/'  ENTER  2  IF  OPTICS  DIFFRACTION  IS  TO  BE  INPUTTED') 

READ(*,*)INOP 

WRITE(*, 1110) (IHEAD(I), 1=41, 80) 
1110  FORMAT(/'  ENTER  OPTICS  INFORMATION  -  80  CHARACTERS »/2X,40A2,/) 

READ(*, 1030) (IHEAD(I), 1=41, 80) 

11  IF  (INPT-2)  12,20,65 

C  CALCULATION  OF  SOURCE  PATTERN  IF  NOT  INPUTTED 

12  WRITE(*,1120)A 

1120  FORMAT(/'  ENTER  MAXIMUM  OF  SOURCE  IN  VATTS/STERADIAN  ' ,2X,E10.4,/) 

READ(*,*)A 

WRITE(*  1130)B 
1130  FORMATf/'  ENTER  SIGMA  OF  GAUSSIAN  SOURCE  IN  MICRORADIANS' ,2X,E10. 4 

$,/) 

READ(*,*)B 
VRITE(*,1140) 
1140  FORMAT(/'  ENTER  1  IF  LONG- TERM (EXPOSURE)  IS  TO  BE  TAKEN  FOR  OUTPUT 
$'/'  ENTER  2  IF  SHORT- TERM (EXPOSURE)  IS  TO  BE  TAKEN  FOR  OUTPUT') 
READ(*,*)  EXPO 
VRITE(*,1150) 
1150  FORMAT(/'  ENTER  THE  FRACTION  OF  TOTAL  POWER  THAT  IS  REQUIRED  AT  TH 
$E  TARGET') 
READ(*,*)  FRAP 
WRITE(*,1160) 
1160  FORMAT(/'  ENTER  1  IF  ANY  VARIABLE  IS  TO  BE  VARIED'/'  ENTER  2  IF  AN 
$Y  VARIABLE  IS  NOT  TO  BE  VARIED  ') 
READ (*,*) VARY 
IF(VARY-l)  110,110,100 
110  WRITE(*,1170) 
1170  FORMAT(/'  4  VARIABLES  CAN  BE  VARIED  IN  THIS  PROGRAM'/'  SELECT  THE 
SNUMBER  THAT  IS  TO  BE  VARIED'//'  1  -  SQCN'/'  2  -  EXT'/'  3  -  RANGE'/ 
$'  4  -  OPT(DIAMETER)') 
READ(*,*)SEL 

GOTO  (120,130,140,150)  SEL 
120  WRITE(*,1180) 
1180  FORMAT(/'  ENTER  THE  MIN-VALUE  OF  SQCN') 
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READ(*,*)SQCNL 

VRITE(*,1190) 
1190  F0RMAT(/'  ENTER  THE  MAX- VALUE  OF  SQCN') 

READ(*,*)SQCNH 

CALL  GETVAR(SQCNL,SQCNH,VA) 

INU=1 
200  SQCN=VA(INU) 

GO  TO  100 
130  VRITE(*,1200) 
1200  FORMAT(/'  ENTER  THE  MIN- VALUE  OF  EXT') 

READ(*,*)EXTL 

WRITE(*,1210) 
1210  FORMAT (/'  ENTER  THE  MAX- VALUE  OF  EXT') 

READ(*,*)EXTH 

CALL  GETVAR(EXTL,EXTH,VA) 

INU=1 
300  EXT=VA(INU) 

GO  TO  100  ' 
140  WRITE(*,1220) 
1220  FORMAT(/'  ENTER  THE  MIN- VALUE  OF  RANGE') 

READ(*,*)RANGL 

WRITE (*, 1230) 
1230  FORMAT (/'  ENTER  THE  MAX- VALUE  OF  RANGE') 

READ(*,*)RANGH 

CALL  GETVAR(RANGL,RANGH,VA) 

INU=1 
400  RANGE=VA(INU) 

GO  TO  100 
150  WRITE(*,1240) 
1240  FORMAT(/'  ENTER  THE  MIN- VALUE  OF  OPT(DIAMETER) ') 

READ(*,*)OPTL 

WRITEf*  1250) 
1250  FORMAT(/'  ENTER  THE  MAX- VALUE  OF  OPT(DIAMETER) ') 

READ(*,*)OPTH 

CALL  GETVAR(OPTL,OPTH,VA) 

INU=1 
500  OPT=VA(INU) 
100  IF(VARY.EQ.2)  INU=20 

VRITE(*,1260)  INU,VA(INU) 
1260  FORMAT(/'  CONTINUE' ,2X, 12, 2X,E10. 4) 

DO  15  1=1,512 

F=SCALE*FL0AT(I-1) 

ARG=F*F/(2.*B*B) 

IF(ARG.GT.70)  GO  TO  14 

RG1(I)=A*EXP(-ARG) 

K=I 

GO  TO  15 

14  RG1(I)=0 

15  CONTINUE 

CALL  0UTIN(RG1,1,INU,2) 
GO  TO  25 
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C  INPUT  OF  SOURCE  PATTERN  IF  NOT  CALCULATED 
20  CALL  OUTIN(RG1,0,20,2) 
DO  22  1=1,512 
K=513- 1 

IF(RGl(K).NE.O.)  GO  TO  25 
22  CONTINUE 
C  SOURCE  PSF  IS  CONVERTED  TO  LSF 
25  CALL  LSF(RG1,K) 

CALL  OUTIN(RG1,1,INU,10) 
C  FOURIER  TRANSFORM  OF  SOURCE  LSF 
CALL  FFT(RG1,1,K) 
CALL  0UTIN(RG1,1,INU,18) 
IF(IN0P- 1)30,30,40 
C  CALCULATION  OF  OPTICS  DIFFRACTION  PATTERN  IF  NOT  INPUTTED 
30  OBS=OBRA*OBRA 
0B=(1.-0BS)**2 
Q=3.14159E-6*SCALE*0PT/WAV 
DO  34  1=1,512 
R=Q*FL0AT(I-1) 
IF(R.GT.60.)  GO  TO  32 
RG2(I)=((AIRY(R)-0BS*AIRY(R*0BRA))**2)/0B 
K=I 

GO  TO  34 
32  RG2(I)=0. 
34  CONTINUE 

CALL  0UTIN(RG2,1,INU,26) 
GO  TO  46 
C  INPUT  OF  OPTICS  DIFFRACTION  IF  NOT  CALCULATED 
40  CALL  OUTIN(RG2,0,20,26) 
DO  44  1=1,512 
K=512-I 

IF(RG2(K).NE.O.)  GO  TO  46 
44  CONTINUE 
C  CONVERSION  OF  OPTICS  DIFFRACTION  PATTERN  (PSF)  TO  OPTICS  LINE  SPREAD 
C  LSF 

46  CALL  LSF(RG2,K) 

CALL  0UTIN(RG2,1,INU,34) 
C  FOURIER  TRANSFORM  OF  OPTICS  LSF 

CALL  FFT(RG2,1,K) 
C  NORMALIZATION  OF  THE  FOURIER  TRANSFORM  OF  OPTICS  LSF 
DIV=RG2(1) 
DO  57  1=1,512 
57  RG2(I)=RG2(I)/DIV 

CALL  0UTIN(RG2,1,INU,42) 
C  PRODUCT  OF  FOURIER  TRANSFORM  OF  SOURCE  LSF  WITH  FOURIER  TRANSFORM  OF 
C  LSF 

DO  64  1=1,512 
64  RG1(I)=RG1(I)*RG2(I) 
CALL  OUTIN(RG1,1,INU,50) 
GO  TO  66 
C  INPUT  OF  FOURIER  TRANSFORM  OF  SOURCE  LSF  *  FOURIER  TRANSFORM  OF  OPTIC 
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65  CALL  OUTIN(RG1,0,20,50) 

CALCULATION  OF  ATMOSPHERIC  MODULATION  TRANSFER  FUNCTION  (LONG  TERM) 

66  CALL  AMTF(O) 

CALL  0UTIN(RG2,1,INU,58) 
PRODUCT  OF  AMTF  WITH  FFT  OF  SOURCE  LSF  AND  FFT  OF  OPTICS  LSF  (LONG  T 

DO  68  1=1,512 
68  RG2(I)=RG1(I)*RG2(I) 

CALL  OUTIN(RG2,l,INU,66) 
INVERSE  FOURIER  TRANSFORM  OF  IMAGE  LSF  (LONG  TERM) 

CALL  FFT(RG2,0,K) 

CALL  OUTIN(RG2,l,INU,74) 
CONVERSION  OF  IMAGE  LSF  TO  PSF  (LONG  TERM) 


CALL  ABEL(K) 

J(RG2,1 
C  CALCULATION  OF  FRACTION  OF  POWER  TO  TARGET  INSIDE  A  CIRCLE  OF  RADIUS 


CALL  0UTIN(RG2,1,INU,82) 


C  (LONG  TERM) 

CALL  INT(K,TPOW) 
IF (EXPO. EQ. 2)  GOTO  72 
CALL  RESUL(INU,TPOW,RADIU,POWER,FRAP) 
WRITE(*,1270)  TPOW,RADIU(INU) , FRAP, POWER (INU) 
1270  FORMAT(/'  TOTAL  P=' ,E10.4,2X, »RADIUS=' ,E10.4,2X, 'FRAC=' ,F5.3,2X, ' 
SAVE  P=',E10.4) 
72  CALL  OUTIN(RG2,1,INU,90) 
C  CALCULATION  OF  ATMOSPHERIC  MODULATION  TRANSFER  FUNCTION  (SHORT  TERM) 
CALL  AMTF(l) 

CALL  0UTIN(RG2,1,INU,98) 
C  PRODUCT  OF  AMTF  WITH  FFT  OF  SOURCE  LSF  AND  FFT  OF  OPTICS  LSF  (SHORT 
DO  74  1=1,512 
74  RG2(I)=RG1(I)*RG2(I) 
CALL  OUTIN(RG2,1,INU,106) 
C  INVERSE  FOURIER  TRANSFORM  OF  IMAGE  LSF  (SHORT  TERM) 
CALL  FFT(RG2,0,K) 
CALL  0UTIN(RG2,1,INU,114) 
C  CONVERSION  OF  IMAGE  LSF  TO  PSF  (SHORT  TERM) 
CALL  ABEL(K) 

CALL  0UTIN(RG2,1,INU,122) 
C  CALCULATION  OF  FRACTION  OF  POWER  TO  TARGET  INSIDE  A  CIRCLE  OF  RADIUS 
C  (SHORT  TERM) 

CALL  INT(K,TPOW) 
IF(EXPO.Eq.l)  GOTO  76 
CALL  RESUL(INU,TPOW,RADIU,POWER,FRAP) 
WRITE(*,1270)  TPOW,RADIU(INU),FRAP,POWER(INU) 
76  CALL  OUTIN(RG2,1,INU,130) 
INU=INU+1 
IFfSEL.EQ.l. AND.INU.LT. 21)  THEN 

GOTO  200 
ELSEIF(SEL.EQ. 2. AND.INU.LT. 21)  THEN 

GOTO  300 
ELSEIF(SEL.EQ. 3. AND.INU.LT. 21)  THEN 

GOTO  400 
ELSEIF(SEL.EQ. 4. AND.INU.LT. 21)  THEN 
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GOTO  500 
ELSE 

GOTO  80 
ENDIF 
80  CALL  VRITO(VA,RADIU, POWER, ABC) 
CALL  OUTIN(ABC,2,20,1) 
STOP 
END 


SUBROUTINE  RESUL(INU,TPOW,RADIU,POWER,FRAP) 

DIMENSION  RADIU(2O),P0WER(2O) 

C0MM0N/REAL/RG1(512),RG2(512) 

COMMON//WAV, OPT, OBRA, SCALE, RANGE, SqCN, EXT, A, B,IDATE(30),IHEAD(80) 

DO  5  1=1,512 

P=RG2(I)/FRAP 

IF(P.GE.l.O)  GOTO  7 
5  CONTINUE 
7  RADIU(INU)=RANGE*FL0AT(I)*SCALE*l.E-6 

P0WER(INU)=(TP0W*FRAP)/(3.14159*RADIU(INU)**2) 

RETURN 

END 


SUBROUTINE  VRITO ( VA , RADIU, POWER, ABC) 

INTEGER*2  SEL,VARY 

DIMENSION  VA(20) ,RADIU(20) ,P0WER(2O) ,ABC(64) 

COMMON/OJ/SEL,VARY 

IF(VARY.EQ.2)  GOTO  100 

DO  10  1=1,20 

WRITE(9,20)  VA(I) ,RADIU(I) ,POWER(I) 
20  FORMAT(2X,E10.4,2X,F10.3,2X,E10.4) 
10  CONTINUE 

WRITE(9,30) 
30  FORMAT (2X,1  VARIABLE' ,7X, 'RADIUS' ,3X, 'AVERAGE  USABLE  POWER') 

GOTO  (2,3,4,5)SEL 

2  WRITE(9,40) 

40  F0RMAT(/,2X, 'VARIABLE  IS  SQCN',/) 
GO  TO  200 

3  WRITE(9,50) 

50  F0RMAT(/,2X, 'VARIABLE  IS  EXT',/) 
GO  TO  200 

4  WRITE (9, 60) 

60  F0RMAT(/,2X, 'VARIABLE  IS  RANGE',/) 
GO  TO  200 

5  WRITE(9,70) 

70  FORMAT (/,2X, 'VARIABLE  IS  OPT(DIAMETER) ' ,/) 

GO  TO  200 
100  WRITE(9,80)  RADIU(20) ,P0WER(2O) 
80  FORMAT(2X,F10.3,2X,E10.4,/) 

WRITE(9,90) 
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90  F0RMAT(2X, 'RADIUS' ,6X, 'AVERAGE  USABLE  POWER') 
200  WRITE(9,300)(ABC(I),I=1,64) 
300  FORMAT ('WAVELENGTH  =  '  E12.5  »  DIAMETER  OF  OPTICS  =  * 


SE12.5,/, 'OBSCURATION  RATIO     =  ',E12.5,'  SCALE  =  ' 

$,E12.5,/,*RANGE  =  ',E12.5,'  CN  SQUARED 

$',E12. 5, //EXTINCTION  COEEFICIENT  =  ',E12.5,*  AMPLITUDE 

$  ',E12.5,/,'SIGMA  =  »,E12.5,/,2X,15A4,/,2X,20A4,/ 

$,2X,20A4) 

RETURN 

END 


SUBROUTINE  GETVAR(MIN,MAX,VA) 
DIMENSION  VA(20) 
REAL  MIN,MAX,K,L 
K=MAX-MIN 
L=K/19. 
VA(1)=MIN 
DO  10  1=2,19 
VA(I)=VA(1)+L*FL0AT(M) 
10  CONTINUE 
VA(20)=MAX 
RETURN 
END 


SUBROUTINE  ABEL  (K) 
CONVERTS  LINE  SPREAD  FUNCTION  (X)  TO  POINT  SPREAD  FUNCTION  (R) 
C0MM0N/REAL/RG1(512),RG2(512) 
G=RG2(1) 
F=G 

RG2(1)=.2*(RG2(1)-RG2(2)) 
DO  1  1=1,510 
H=RG2(I+1) 

RG2(I+1)=.4*G+.2*H-  .6*RG2(I+2) 
C=\\ 
IF(ABS(G).GT.1E-6*F)  K=I 

1  CONTINUE 
DO  3  1=1, K 

RG2(I)=RG2(I)/(2.*SQRT((FL0AT(I)-.9)**2-(FL0AT(I-l))**2)) 
11=1+1 

DO  2  J=II,K 
RG2(I)=RG2(I)+RG2(J)/SQRT((FL0AT(J)-.9)**2-(FL0AT(I-l))**2) 

2  CONTINUE 
RG2(I)=RG2(I)/3. 14159 

3  CONTINUE 
K=K+1 

DO  4  I=K,512 

4  RG2(I)=0. 
RETURN 
END 
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i)  GO  TO  1 

<X-X*X/3 


FUNCTION  AIRY(X) 
AIRY(X)  =  2*J1(X)/X  WHERE  J1(X)  IS  THE  BESSEL  FUNCTION  OF  ORDER  ONE 
IF(X.EQ.O.)  GO  TO  8 
2=0 

IF(X.GT.15' 
J=20.+10.*J 
GO  TO  2 

1  J=90.+X/2. 

2  IF(X.LT.5.)  GO  TO  3 
W=1.4*X+60./X 
GO  TO  4 

3  W=X+6. 

4  I=MAXl(V,(3.+X/4.)) 
DO  6  M=I,J,3 
V=l.E-28 
S=0. 
T=S 

F=l.-4.*(FL0AT(M)/2.-FL0AT(M/2)) 
MM=M-2 
DO  5  K=1,MM 
U=2.*FLOAT(M-K)*V/X-T 
T=V 
V=U 
IF((M-K-2).EQ.O)  B=U 

FF=1.-F 

5  S=S+V*FF 
U=2.*V/X-T 

s=s+u 

B=B/S 

IF((ABS(B-Z)-ABS(B*l.E-6)).LT.O.)  GO  TO  7 

6  CONTINUE 

7  AIRY=2.*B/X 
RETURN 

8  AIRY=1. 
RETURN 
END 


SUBROUTINE  AMTF(J) 
C  THIS  SUBROUTINE  CALCULATES  THE  ATMOSPHERIC  MODULATION  TRANSFER  FUNCTI 
C  FOR  LONG  TERM  (J=0)  AND  SHORT  TERM  (J=l)  AND  STORES  THE  RESULT  IN  RG2 
C  THE  RESULT  IS  THEN  MULTIPLIED  BY  EXTINCTION  FACTOR  EXP(-EXT*RANGE) . 

C0MM0N/REAL/RG1(512),RG2(512) 

COMMON//WAV, OPT, OBRA, SCALE, RANGE, SQCN, EXT, A, B,IDATE(30) ,IHEAD(80) 

AA=21 . 6*SQCN*RANGE*WAV** (- . 3333333) 

EX=EXP(-EXT*RANGE) 

D=0. 

RG2(1)=EX 

L=0 
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DO  20  1=2,512 

IF(L)  10,5,15 
5  FF=976.56*(M)/SCALE 

ARG=AA*(FF**1. 666667) 

IF(J.EQ.l)  ARG=ARG*(1-  (FF*WAV/OPT)**. 3333333) 

IF(ARG.GT.13)  L=l 

IF((ARG-D).LT.O.)  L=- 1 

RG2(I)=EX*EXP(-ARG) 

D=ARG 

M=I 

GO  TO  20 
10  RG2(I)=RG2(M)*(C0S((I-M)*3.14159/M))**2 

IF(I.LT.3*M/2)  GO  TO  20 

L=l 
15  RG2(I)=0. 
20  CONTINUE 

RETURN 

END 


SUBROUTINE  BITINV(K,L,M) 
C  FINDS  THE  BIT  INVERSION  OF  A  (M)  BIT  NUMBER  (J)  TO  (L) 
C  THIS  IS  USED  IN  THE  FOURIER  TRANSFORM  SUBROUTINE 

J=K 

L=0 

DO  1  1=1, M 

L=L*2 

IF(J.NE.J/2*2)  L=L+1 
1  J=J/2 

RETURN 

END 


SUBROUTINE  CL0SE(IE1,IE2,IE3) 

INTEGER*2  LU,LI,IH,IM,IS,IHS,IY,IMO,ID,SEL,VARY 

COMMON/OI/LU,LI,IH,IM,IS,IHS,IY,IMO,ID 

CL0SE(LU,I0STAT=IE1 

CL0SE(LI,I0STAT=IE2( 

CL0SE(9,I0STAT=IE3) 

RETURN 

END 


SUBROUTINE  VRITF(IER,BUF,L) 

DIMENSION  BUF(512) 

INTEGER*2  LU,LI,IH,IM,IS,IHS,IY,IMO,ID,SEL,VARY 

COMMON/OI/LU,LI,IH,IM,IS,IHS,IY,IMO,ID 

1  FORMAT(lHl) 

2  FORMATflX,/, IX, //"GAUSSIAN  SOURCE  FUNCTION'") 

3  F0RMAT(2X,5E12.4) 

4  F0RMAT(5X,I5,/,1X) 
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10  F0RMAT(1X,/,1X,/, * ' 'SOURCE  LINE  SPREAD  FUNCTION"') 

18  F0RMAT(1X,/,1X,/,'"FFT  OF  LINE  SPREAD  FUNCTION  OF  SOURCE'") 

26  F0RMAT(1X,/,1X,/,"  'OPTICS  DIFFRACTION  PATTERN  0(R)»") 
34  F0RMAT(1X,/,1X,/,'"0PTIS  LINE  SPREAD  FUNCTION  OLSF(X)'" 
42  F0RMAT(1X,/,1X,/,'"FFT  OF  OPTICS  LINE  SPREAD  FUNCTION'" 
50  F0RMAT(1X,/,1X,/,,MFFT  OF  SOURCE  *  FFT  OF  OPTICS'") 
58  F0RMAT(1X,/,1X,/,'"L0NG  TERM  MTF  OF  THE  ATMOSPHERE"') 

66  F0RMAT(1X/,1X/, ' ' 'FFT  OF  IMAGE  LINE  SPREAD  FUNCTION  LONG  TERM'  ") 
74  F0RMAT(1X, /, IX, /,' "IMAGE  LINE  SPREAD  FUNCTION  LONG  TERM'") 
82  F0RMAT(1X, /, IX, /,' "IMAGE  POINT  SPREAD  FUNCTION  LONG  TERM'") 
90  FORMAT  MX/,  IX, /,"  'ENERGY  INSIDE  CIRCLE  OF  RADIUS  R  LONG  TERM*") 
98  F0RMAT(1X, /, IX, /,' "SHORT  TERM  ATMOSPHERIC  MTF'") 
106  F0RMAT(1X/1X/, ' * 'FFT  OF  IMAGE  LINE  SPREAD  FUNCTION  SHORT  TERM'  ") 
114  FORMAT ( IX, /, IX, /,' "IMAGE  LINE  SPREAD  FUNCTION  SHORT  TERM'") 
122  FORMATMX, /, IX, /,' "IMAGE  POINT  SPREAD  FUNCTION  SHORT  TERM'") 
130  F0RMAT(1X/,1X/, ' ' 'ENERGY  INSIDE  CIRCLE  OF  RADIUS  R  SHORT  TERM'  ") 
150  F0RMAT(2X,2I2, ' : ' ,12, ' : ' ,I2,2X,I2, '- ' ,12, '- ' ,14) 
IER=0 

WRITE(LU,1) 
IF(L.EQ.l)  GO  TO  20 
WRITE(LU,4)  L 
WRITE(LU,3)  BUF 
J=(L+6)/8  ' 

GO  TO  (5, 7, 9, 11, 13, 15, 17, 19, 21, 23, 25, 27, 29, 31, 33, 35, 37), J 
5  WRITE(LU,2) 

GO  TO  38 
7  VRITE(LU,10) 

GO  TO  38 
9  WRITE(LU,18) 
GO  TO  38 

11  WRITE(LU,26) 
GO  TO  38 

13  VRITE(LU,34) 

GO  TO  38 
15  WRITE(LU,42) 

GO  TO  38 
17  WRITE(LU,50) 

GO  TO  38 

19  WRITE(LU,58) 
GO  TO  38 

21  WRITE(LU,66) 

GO  TO  38 
23  WRITE (LU, 74) 

GO  TO  38 
25  WRITE(LU,82) 

GO  TO  38 

27  VRITE(LU,90) 
GO  TO  38 

29  WRITE(LU,98) 

GO  TO  38 
31  WRITE(LU,106) 

43 


GO  TO  38 
33  VRITE(LU,114) 

GO  TO  38 
35  WRITE(LU,122) 

GO  TO  38 

37  WRITE(LU,130) 

38  WRITE(LU,150)IH,IM,IS,IHS,IMO,ID,IY 
RETURN 

20  VRITE(LU,6)(BUF(I),I=1,64) 
6  FORMAT ('WAVELENGTH  =  »,E12.5,'  DIAMETER  OF  OPTICS  =  ' 


$E12. 5, //OBSCURATION  RATIO     =  '  E12.5,'  SCALE  =  ' 

$,E12.5,/,'RANGE  =  ',E12.5,'  CN  SQUARED 

$»,E12. 5, //EXTINCTION  COEEFICIENT  =  »,E12.5,'   AMPLITUDE 

$  ' ,E12.5,/,'SIGMA  =  ,,E12.5,/,2X,15A4,/,2X,20A4,/ 

$,2X,20A4) 

RETURN 

END 


SUBROUTINE  FFT(RG,ID,KA) 
FOURIER  TRANSFORM  SUBROUTINE 

IMPLICIT  COMPLEX  (C) 

DIMENSION  CE(512),CG(1024) 

DIMENSION  RG(512) 

DATA  1ST, TH/O,. 006135923/ 

CG(1)=RG(1) 

DO  10  1=2,512 

CG(I)=RG(I) 
10  CG(1026-I)=CG(I) 

CG(513)=CG(512) 

IF(IST.NE.O)  GO  TO  2 

IST=1 

DO  1  1=1,512 

A=TH*FL0AT(I-1) 

C=CMPLX(0.,A) 

1  CE(I)=CEXP(C) 
CE(257)=(0.,1.) 

2  DO  5  M=l,10 
K=2**(10-M) 
L=2**(M-1) 
DO  4  LL=1,L 
KC=LL- 1 

CALL  BITINV  (KC,N,9) 

CD=CE(N+1) 

IF(ID.EQ.l)  CD=CONJG(CD) 

KL=2*KC*K+1 

KU=KL+K- 1 

DO  3  I=KL,KU 

CA=CG(I)+CG(I+K)*CD 

CG(I+K)=CG(I)-CG(I+K)*CD 

3  CG(I)=CA 
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4  CONTINUE 

5  CONTINUE 

DO  8  1=1,1024 

J=I-1 

CALL  BITINV  (J,N,10) 

IF(J-N)  6,7,8 

6  N=N+1 
CA=CG(I)/(32.,0.) 
CG(I)=CG(N) 
CG(N)=CA 

7  CG(I)=CG(I)/(32.,0.) 

8  CONTINUE 

DO  20  1=1,512 
RG(I)=CG(I) 

IF(RG(I).GT.1.E-6*RG(1))  KA=I 
20  CONTINUE 
RETURN 
END 


SUBROUTINE  INT(K,TPOW) 
C  CALCULATES  THE  FRACTION  OF  POWER  INSIDE  A  CIRCLE  OF  RADIUS  R 
C  OR  CALCULATES  THE  POWER  INSIDE  A  CIRCLE  OF  RADIUS  R 

C0MM0N/REAL/RG1(512),RG2(512) 

COMMON/ /WAV , OPT , OBRA , SCALE , RANGE , SQCN , EXT , A , B , IDATE (30) , IHEAD (80) 

SR=(SCALE*RANGE)**2 

RG2(1)=.7853982E- 12*RG2(1)*SR 

DO  3  1=2,512 

IF(I-K)  1,1,2 

1  RG2(I)=RG2(I-1)+6.28315E-12*RG2(I)*SR*FL0AT(I-1) 
GO  TO  3 

2  RG2(I)=RG2(K) 

3  CONTINUE 
TP0W=RG2(512) 
DO  5  1=1,512 

5  RG2(I)=RG2(I)/TP0W 
RETURN 
END 


SUBROUTINE  LSF  (RG3,K) 
CONVERTS  POINT  SPREAD  FUNCTION  (R)  TO  LINE  SPREAD  FUNCTION  (X) 
DIMENSION  RG3(512) 
DO  2  1=1, K 
Z=RG3(I) 
J=l 
1  R=SQRT(FL0AT((I-1)**2+J*J)) 
IR=IFIX(R) 
RI=FLOAT(IR) 
RF=R-  RI 
Z=Z+2*((l.-RF)*RG3(IR+l)+RF*RG3(IR+2)) 
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J=J+1 

IF  (IR.LT.K)  GO  TO  1 

RG3(I)=Z 

CONTINUE 

RETURN 

END 


SUBROUTINE  OPEN 

INTEGER*2  LU,LI,IH,IM,IS,IHS,IY,IMO,ID,SEL,VARY 

COMMON/OI/LU,LI,IH,IM,IS,IHS,IY,IMO,ID 

OPEN(7,IOSTAT=IE1,ERR=10,FILE='B:INPUT.DTA',STATUS='OLD') 

OPEN (8, I0STAT=IE2,ERR=11,FILE='B: OUTPUT. DTA') 

0PEN(9,I0STAT=IE3,ERR=12,FILE=,B:0Ur2.DTA') 

RETURN 

10  VRITE(*,13)IE1 

13  FORMAT ('  ERROR  STOP  IN  OPENING  INPUT  FILE  ERN  =  ',12) 
STOP 

11  VRITE(*,14)IE2 

14  FORMAT (■  ERROR  STOP  IN  OPENING  0UTPUT1  FILE  ERN  =  ',12) 
STOP 

12  VRITE(*  15)IE3 

15  FORMAT^  ERROR  STOP  IN  OPENING  0UTPUT2  FILE  ERN  =  ',12) 
END 


SUBROUTINE  OUTIN(BUF,J,K,L) 

INTEGER*2  NAME,LU,LI,IH,IM,IS,IHS,IY,IMO,ID,SEL,VARY 

COMMON/OI/LU,LI,IH,IM,IS,IHS,IY,IMO,ID 

DIMENSION  BUF(1024) 

CALL  GETTIM  (IH,IM,IS,IHS) 

IF(J)  20,40,60 
20  CALL  OPEN 

RETURN 
40  CALL  READF(IER,BUF,L) 

IF(IER.GE.O)  RETURN 

WRITE(LU,1005)IER,IN 
1005  FORMAT (■  ERROR  STOP.  READF  ERROR.  IER  =  ',16, 'BLOCK  #  =  ',16) 

STOP 
60  IF(K.LT.20)RETURN 

CALL  VRITF(IER,BUF,L) 

IF(J.EQ.2.AND.IER.GE.O)  GO  TO  70 

IF(IER.GE.O)  RETURN 

WRITE(LU,1140)IER,IN 
1140  FORMAT('  ERROR  STOP.  WRITF  ERROR.  IER  =  ',16,'  BLOCK  #  =  ',16) 

STOP 
70  VRITE(*,1150) 

CALL  CL0SE(IE1,IE2,IE3) 
1150  FORMAT ('  CALCULATION  COMPLETE') 

STOP 

END 
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SUBROUTINE  READF(IER,BUF,L) 

INTEGER*2  LU,LI,IH,IM,IS,IHS,IY,IMO,ID,SEL,VARY 

DIMENSION  BUF(512) 

COMMON/OI/LU,LI,IH,IM,IS,IHS,IY,IMO,ID 

IF(L.EQ.l)  GO  TO  20 

1  FORMAT(8F10.5) 
DO  10  1=1,512,8 
K=I+7 
READ(LI,1)(BUF(J),J=I,K) 

10  CONTINUE 

IER=1 

RETURN 
20  READ(LI,2)(BUF(I),I=1,9) 

2  FORMAT (3F10. 5) 
READ(LI,3)(BUF(I),I=10,24) 

3  F0RMAT(15A4) 
READ(LI,4)(BUF(I),I=25,64) 

4  FORMAT(20A4) 
RETURN 

END 
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APPENDIX  B 


A.   MODULE  INTERCONNECTIONS 


GETVAR 

AIRY 

LSF 

AMTF 

FFT 

BITINV 

MAIN 

ABEL 

OPEN 

INT 

READF 

OUTIN 

WRITF 

WRITO 

CLOSE 

RESUL 

[Fig.  4]      Module  interconnections 
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B.      INPUT  FILE 


1.06e-06  0.095  0.000 

1.00000        2000.0    5.00e-14 
0.00005      50000.0        1.5000 
H0UR:MIN  DATE  MONTH  YEAR 
SECOND  LINE  OF  COMMENTS 
THIRD  LINE  OF  COMMENTS 


1.0 

1.0 

1.0 

1.0 

1.0 

1.0 

1.0 

1.0 

1.0 

1.0 

1.0 

0.8 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0,0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 
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0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 
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APPENDIX  C 
OUTPUT  OF  THE  SAMPLE  RUNS 


OUTPUT  OF  RUN  1 


VARIABLE 

RADIUS 

AVERAGE  U 

.8885E-14 

.240 

.9770E+00 

.2030E-13 

.360 

.4325E+00 

.3172E-13 

.470 

.2527E+00 

.4313E-13 

.550 

. 1838E+00 

.5455E-13 

.630 

.1395E+00 

.6597E-13 

.700 

.1125E+00 

.7738E-13 

.770 

.9253E-01 

.8880E-13 

.830 

.7925E-01 

.1002E-12 

.880 

.7015E-01 

.1116E-12 

.940 

.6115E-01 

.  1231E- 12 

.990 

.5483E-01 

.  1345E- 12 

1.040 

.4939E-01 

.1459E-12 

1.080 

.4553E-01 

.1573E-12 

1.130 

.4132E-01 

.1687E-12 

1.170 

.3829E-01 

.1801E-12 

1.210 

.3555E-01 

.1916E-12 

1.250 

.3306E-01 

.2030E-12 

1.290 

.3080E-01 

.2144E-12 

1.320 

.2918E-01 

.2258E-12 

1.360 

.2726E-01 

VARIABLE  IS  SQCN 

WAVELENGTH 

OBSCURATION  RATIO 

RANGE 

EXTINCTION  COEEFICIENT  = 

SIGMA 
05:18  19  12  1989 

Assumed  Gaussian  source  for  simulation 
Assumed  optics  information  for  run  1 


10600E-05  DIAMETER  OF  OPTICS  =  .95000E-01 

00000E+00  SCALE  =  .10000E+01 

10000E+05  CN  SQUARED  =  varied 

50000E-04  AMPLITUDE  =  .50000E+03 
15000E+01 
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B.      OUTPUT  OF  RUN  2 


VARIABLE 


RADIUS      AVERAGE  USABLE  POWER 


.5000E-05 
.3105E-04 
.5711E-04 
.8316E-04 
.1092E-03 
.1353E-03 
.1613E-03 
.1874E-03 
.2134E-03 
.2395E-03 
.2655E-03 
.2916E-03 
.3176E-03 
.3437E-03 
.3697E-03 
.3958E-03 
.4218E-03 
.4479E-03 
.4739E-03 
.5000E-03 

VARIABLE  IS  EXT 

WAVELENGTH 

OBSCURATION  RATIO 

RANGE 

EXTINCTION  COEEFICIENT  = 

SIGMA 
06:15  19  12  1989 

Assumed  Gaussian  source  for  simulation 
Assumed  optics  information  for  run  2 


.250 

.1411E+01 

.250 

.1088E+01 

.250 

.8383E+00 

.250 

.6460E+00 

.250 

.4978E+00 

.250 

.3837E+00 

.250 

.2957E+00 

.250 

.2279E+00 

.250 

.1756E+00 

.250 

. 1353E+00 

.250 

.1043E+00 

.250 

.8037E-01 

.250 

.6194E-01 

.250 

.4773E-01 

.250 

.3678E-01 

.250 

.2835E-01 

.250 

.2185E-01 

.250 

.1683E-01 

.250 

.1297E-01 

.250 

.9998E-02 

.10600E-05  DIAMETER  OF  OPTICS  = 

.95000E-01 

= 

.OOOOOE+OO  SCALE 

.10000E+01 

= 

.10000E+05  CN  SQUARED 

.  10000E- 13 

NT  = 

varied   AMPLITUDE 

.50000E+03 

= 

.15000E+01 
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C.   OUTPUT  OF  RUN  3 


VARIABLE 


RADIUS   AVERAGE  USABLE  POWER 


.5000E+04 
.5526E+04 
.6053E+04 
.6579E+04 
.7105E+04 
.7632E+04 
.8158E+04 
.8684E+04 
.9211E+04 
.9737E+04 
. 1026E+05 
.1079E+05 
.1132E+05 
.1184E+05 
.1237E+05 
.1289E+05 
.1342E+05 
.1395E+05 
. 1447E+05 
.1500E+05 

VARIABLE  IS  RANGE 

WAVELENGTH 

OBSCURATION  RATIO 

RANGE 

EXTINCTION  COEEFICIENT  = 

SIGMA 

10:58  19  12  1989 

Assumed  Gaussian  source  for  simulation 

Assumed  optics  information  for  run  3 


090 

.2234E+01 

105 

.1952E+01 

115 

.1901E+01 

132 

.1671E+01 

149 

.1476E+01 

168 

. 1309E+01 

188 

.1167E+01 

200 

.1136E+01 

221 

.1016E+01 

243 

.9120E+00 

257 

.8882E+00 

281 

.7997E+00 

306 

.7221E+00 

320 

.7033E+00 

346 

.6368E+00 

374 

.5781E+00 

389 

.5630E+00 

418 

.5124E+00 

434 

.4990E+00 

465 

.4551E+00 

.10600E-05 

DIAMETER  OF  OPTICS  = 

.95000E-01 

= 

.00000E+00 

SCALE 

.10000E+01 

= 

varied 

CN  SQUARED 

.10000E-13 

T  = 

.50000E-04 

AMPLITUDE 

.50000E+03 

= 

.15000E+01 
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D.   OUTPUT  OF  RUN  4 


VARIABLE 

RADIUS 

AVERAGE  U 

.2500E-01 

.330 

.5145E+00 

.3026E-01 

.310 

.5840E+00 

.3553E-01 

.290 

.6672E+00 

.4079E-01 

.280 

.7166E+00 

.4605E-01 

.280 

.7173E+00 

.5132E-01 

.270 

.7715E+00 

.5658E-01 

.270 

.7716E+00 

.6184E-01 

.260 

.8321E+00 

.6711E-01 

.260 

.8321E+00 

.7237E-01 

.260 

.8321E+00 

.7763E-01 

.260 

.8321E+00 

.8289B-01 

.250 

.9000E+00 

.8816E-01 

.250 

.9000E+00 

.9342E-01 

.250 

.9000E+00 

.9868E-01 

.250 

.9000E+00 

. 1039E+00 

.250 

.9000E+00 

. 1092E+00 

.250 

.9000E+00 

.1145E+00 

.250 

.9000E+00 

.1197E+00 

.250 

.9000E+00 

.1250E+00 

.250 

.9000E+00 

VARIABLE  IS  OPT(DIAMETER) 

WAVELENGTH 

OBSCURATION  RATIO 

RANGE 

EXTINCTION  COEEFICIENT  = 

SIGMA 

11:38  19  12  1989 

Assumed  Gaussian  source  for  simulation 

Assumed  optics  information  for  run  4 


10600E-05  DIAMETER  OF  OPTICS  =  varied 

00000E+00  SCALE  =  .10000E+01 

10000E+05  CN  SQUARED  =  .10000E-13 

50000E-04  AMPLITUDE  =  .50000E+03 
15000E+01 
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E.      GRAPHICAL  OUTPUT  OF  CALCULATION  STEPS 
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[Fig.   5]   Gaussian  source  function 
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[Fig.   6]   Source  line  spread  function 
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[Fig.   7]   FFT  of  source  LSF 
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[Fig.   8]   Optics  diffraction  pattern 
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[Fig.  9]   Optics  line  spread  function 
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[Fig.   10]   FFT  of  optics  LSF 


60 


5>  1 

O"    I 


05 


spnjiidtuv 


[Fig.    11]   FFT  of  source  *   FFT  of  optics 
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[Fig.    12]   AMTF  (long  term) 
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[Fig.   13]   FFT  of  image  LSF  (long  term) 
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[Fig.   14]   Image  LSF  (long  term) 
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[Fig.   15]   Image  PSF  (long  term) 
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[Fig.   16]   Energy  inside  circle  of  radius  r  (long  term) 
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[Fig.   17]   AMTF  (short  term) 
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[Fig.    18]   FFT  of  image  LSF  (short   term) 
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[Fig.   19]   Image  LSF  (short  term) 
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[Fig.  20]   Image  PSF  (short  term) 
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[Fig.   21]   Energy  inside  circle  of  radius  r   (short  term) 
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